# This file is part of Sonedyan.
#
# Sonedyan is free software; you can redistribute it and/or
# modify it under the terms of the GNU General Public
# License as published by the Free Software Foundation;
# either version 3 of the License, or (at your option) any
# later version.
#
# Sonedyan is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied
# warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
# PURPOSE.  See the GNU General Public License for more
# details.
#
# You should have received a copy of the GNU General Public.
# If not, see <http://www.gnu.org/licenses/>.
#
# Copyright (C) 2009-2012 Jimmy Dubuisson <jimmy.dubuisson@gmail.com>

# rewire one edge of the graph
#
# mode can be "in" (keep in degrees the same), "out" or "all"
#
# NB: keep the graph strongly connected and forbid self loops
rewire <- function(x, mode = "all")
{
	l <- nrow(x)
	OK <- FALSE
	target <- FALSE
	
	# switch only if old target node has in-degree > 1
	while (!OK)
	{
		vi <- sample(l, 1)
		
		change.in.degrees <- TRUE
		
		if (mode != "out")
		{
			if (mode == "all")
			{
				change.in.degrees <- sample(c(TRUE, FALSE), 1)
			}
			else if (mode == "in")
			{
				change.in.degrees <- FALSE
			}
		}
		
		if (change.in.degrees)
		{
			# let's change an edge target
			nis <- which(x[vi,] > 0)
			vj <- sample(nis, 1)
			
			# we don't want to disconnect the graph
			# in-degree of vj should be > 1
			if (sum(x[,vj]) > 1)
			{
				# choose a new target for vi
				nnis <- which(x[vi,] == 0)
				nnis <- setdiff(nnis, vi)
				lnnis <- length(nnis)
				
				if (lnnis > 0)
				{
					# self-loops are forbidden
					if (lnnis == 1)
					{
						vj2 <- nnis[[1]]
					}
					else
					{
						vj2 <- sample(nnis, 1)
					}
					
					OK <- TRUE
					target <- TRUE
				}
			}
		}
		else
		{
			# let's change an edge source
			nis <- which(x[,vi] > 0)
			vj <- sample(nis, 1)
			
			# we don't want to disconnect the graph
			# out-degree of vj should be > 1
			if (sum(x[vj,]) > 1)
			{
				# choose a new source for vi
				nnis <- which(x[,vi] == 0)
				nnis <- setdiff(nnis, vi)
				lnnis <- length(nnis)
				
				if (lnnis > 0)
				{
					# self-loops are forbidden
					if (lnnis == 1)
					{
						vj2 <- nnis[[1]]
					}
					else
					{
						vj2 <- sample(nnis, 1)
					}
					
					OK <- TRUE
				}
			}
		}
	}
	
	if (target)
	{
		tmp <- x[vi, vj]
		x[vi, vj] <- 0
		x[vi, vj2] <- tmp
	}
	else
	{
		tmp <- x[vj, vi]
		x[vj, vi] <- 0
		x[vj2, vi] <- tmp
	}
	
	return(x)
}

# SANN exponential cooling schedule
exp.cooling.schedule <- function(t, t0 = 1000, alpha = 0.2)
{
	return (t0*alpha^t)
}

# SANN acceptance probability
acceptance.probability <- function(e, enew, T, mode = "min")
{
	if (mode == "min")
	{
		if (enew < e)
		{
			return(1)
		}
		else
		{
			return(exp((e - enew) / T))
		}
	}
	else if (mode == "max")
	{
		if (enew > e)
		{
			return(1)
		}
		else
		{
			return(exp((enew - e) / T))
		}
	}
}

# fn: evaluation function
# gr: new candidates function
# kmax: max number of iterations
# mode: "min" searches for a minimum
optim.sann <- function(mat, fn, gr, kmax, mode = "min")
{
	# Initial state, energy
	s <- mat
	e <- fn(s)
	# initial best solution
	sbest <- s
	ebest <- e
	# energy evaluation count
	k <- 0
	
	# while time left & not good enough
	while (k < kmax)
	{
		# compute new temperature
		T <- exp.cooling.schedule(k/kmax)
		loginfo(paste("Current T°: ", T))
		
		# pick some neighbour
		snew <- gr(s)
		# compute its energy
		enew <- fn(snew)
		loginfo(paste("Current energy: ", enew))
		
		# should we move to it?
		if (acceptance.probability(e, enew, T, mode) > runif(1, 0, 1))
		{
			# yes, change state
			s <- snew
			e <- enew
		}
		
		# is this a new best?
		if ((mode == "min" && e < ebest) || (mode == "max" && e > ebest))
		{
			sbest <- snew
			ebest <- enew
		}
		# increment counter
		k <- k + 1
	}
	
	return(list(s = sbest, e = ebest))
}